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ABSTRACT 


A  convex  nonlinear  programming  problem  (NLP)  can  be 
reformulated  using  an  exact  penalty  function.  One  such 
function  is  the  LI  exact  penalty  function  examined  by  Zangwill 
[11].  The  result  of  this  reformulation  is  an  unconstrained 
optimization  problem.  Unfortunately,  this  nonlinear,  convex 
function  is  not  a  continuously  differentiable  function.  Most 
computer  methods/algorithms  are  based  on  the  condition  that 
the  functions  involved  are  differentiable.  One  exception  is 
the  Ellipsoid  Algorithm,  EA3 ,  by  Ecker  and  Kupferschmid  [9]. 
This  algorithm,  while  demonstrating  only  linear  convergence, 
has  been  shown  to  be  robust  in  solving  problems  with 
nondif ferentiable  functions.  The  problem  reduces  to  the 
following  question.  Is  it  "better"  to  solve  the  original, 
constrained  problem,  or  the  reformulated  unconstrained 
problem?  Here,  "better"  is  defined  as  faster  and  more 
accurate.  After  demonstrating  that  the  exact  penalty 
formulation  will  always  solve  the  convex  NLP,  comparisons  are 
made,  using  the  ellipsoid  algorithm,  of  accuracy  of  solution 
and  solution  time  for  several  example  problems. 
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CHAPTER  l!  THE  EXACT  PENALTY  FUNCTION 


The  convex  nonlinear  programming  problem  seeks  to 
minimize  an  objective  function,  f(x),  for  the  n-dimensional 
vector  y,  subject  to  a  number  of  constraints,  gj  (x)  for 
i=l,2,...,m.  Written  in  standard  form,  the  problem  appears 
as: 

minimize  f(x) 
s.t.  g^(x)  £.0  ,  i-l,2,...,m 

where  f(x)  and  each  of  the  gj(x)  are  convex  functions. 
This  constrained  nonlinear  programming  problem  can  usually  be 
solved  by  one  of  several  algorithms.  This  process,  however, 
can  be  difficult  and  time  consuming.  An  apparently  simpler 
problem  would  be  to  minimize  an  objective  function,  subject  to 
no  constraints.  A  transformation  of  the  constrained  problem 
can  be  made,  resulting  in  an  unconstrained  optimization.  One 
such  reformulation  is  the  exact  penalty  function  presented  by 
Zangwill  in  [11].  The  NLP  reduces  to: 

min  p{x,  c)  ,  xeR  ^ 

where  the  following  functions  and  terms  are  defined: 

pix,  c)  -  fix)  +cPix) 

fix)  ■  original  objective  function 
c  ■  penalty  multiplier 

Fix)  - max  (x)  ,0]  ,  the  exact  penalty  function 
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Penalty  Function  Theory 


The  penalty  function  formulation  introduces  an  increase 
in  objective  function  value  for  infeasible  points,  the 
"penalty”  added  to  the  minimization  problem.  So,  given  a 
particular  solution  point,  x^,  if  a  constraint  is  satisfied, 
gi(x^)<=0,  then  the  max[gj  (x,j) ,  0]=0,  and  the  resulting  penalty 
to  the  objective  value  is  (c)*(0)=0.  On  the  other  hand,  if 
the  constraint  is  violated,  gj(x,^)>0,  then  a  positive 
component  multiplied  by  c  is  added  to  the  objective  function. 

The  constrained  problem  has  now  been  reduced  to  a  single 
unconstrained  minimization.  However,  this  single  function  has 
become  complicated  and  nondif ferentiable.  A  method  will  be 
presented  to  solve  this  reformulated  problem.  First,  it  must 
be  shown  that  this  unconstrained  problem  can  be  solved,  and 
the  solution  of  the  problem  is  the  same  solution  as  for  the 
unconstrained  problem.  Two  assumptions  are  required  for  the 
theorem: 

1.  X*  is  the  optimal  point  for  the  constrained  nonlj.near 
programming  problem. 

2.  5°  €  0  ,  where  S®  is  defined  as  the  interior  of  the 

feasible  set. 

The  penalty  function  theorem  can  then  be  stated  as  follows: 

THEOREM:  Let  X*  be  the  optimal  point  for  a  convex  nonlinear 

programming  problem  and  the  interior  of  the  feasible  set  be 
nonempty.  Then,  the  convex  problem  can  be  solved  by  a  single 
unconstrained  minimization. 

The  proof  of  this  theorem  can  be  found  at  appendix  1.  In  this 
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proof,  the  penalty  multiplier  is  defined  so  that  convergence 
to  the  optimal  solution  is  guaranteed.  The  multiplier  is 
defined  as  follows: 

c  -  P~--  ,  where: 
a 

z  €  S°,  {z  is  a  strictly  feasible  point)  ; 

P  -  fix*)  -fiz)  a  -  max  j  [g^  (z)  ] 

Since  x*,  the  optimal  point,  is  the  solution  being 
sought,  an  estimate  for  f(x*)  must  be  made  in  the  calculation 
of  the  penalty  multiplier. 


Estimation  of  the  Penalty  Multiplier 

The  reformulated  unconstrained  problem  is: 

min  f  (x) +c  ,0] 

However,  c  cannot  be  calculated  directly.  The  following 
method  can  be  used  to  make  an  estimate  of  the  penalty 
multiplier: 

step  1:  find  a  z,  s.t.  gj(z)<0,  ie:  z  is  strictly  feasible. 
step  2 :  calculate  a  =  maXj[gj(z)] 

step  3 :  estimate  P  with: 

P'-f'-f  (z)  ,  where  f'  =  a  lower  bound  on  the  true  objective  value 
step  4 :  the  estimate  for  c  becomes: 

C-  -6^ 

a 
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Now,  any  lower  bound  on  the  true  objective  value  will  yield 
a  c  larger  than  the  true  value  of  c.  The  proof  of  the  penalty 

function  theorem  in  appendix  1  demonstrates  that  any  value 
larger  than  the  true  value  of  c  will  also  cause  the  penalty 
function  formulation  to  converge  to  the  optimal  point. 
Therefore,  the  estimate  for  c  as  calculated  above,  will  be 
sufficiently  large  to  cause  the  problem  to  converge  to  a 
solution. 


An  Illustrative  Example 

The  exact  penalty  function  formulation  and  solution  will 
be  illustrated  by  the  following  example 

min  (Xj^-4)  2+ (>:2-4) - 


s.  t.  Xi-2^0 

X2  -  2  s  0 

-Xj^,  -X2^0 

Graphically,  this  example  can  be  depicted  as  shown  in  figure  1 
on  the  next  page.  The  solution  to  this  simple  convex  problem 
can  be  taken  right  from  the  graph  or  solved  manually  using  the 
KKT  method  described  in  [6].  Obviously,  the  optimal  solution 
to  the  minimization  of  the  radius  of  the  circle  centered  at 
(4,4)  is  x*=(2,2)  with  optimal  value  of  f(x*)=8.  Solving  for 

the  Lagrange  multipliers  yields  the  vector;  >.^-[4  4  0  0] 
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Figure  l 

To  solve  this  problem  using  the  exact  penalty  function 
formulation,  the  problem  will  be  rewritten  as  an  unconstrained 
minimization; 


min  (Xi-4)  2+  (X2-4)  2+c  max  [g^  (x)  ,0] 

and  substituting  each  of  the  four  constraints  into  the  penalty 
function.  The  remaining  problem  is  to  find  a  value  for  the 
penalty  multiplier.  To  begin,  select  the  strictly  feasible 

point,  z=(l,l).  Calculations  yield  f(z)=18  and  a-  -1. 

Substituting  all  of  the  intermediate  terms  into  the  penalty 
multiplier  equation  yields; 

a  -1 

To  complete  the  calculation  of  the  penalty  multiplier,  a  lower 
bound  is  needed  on  the  true  objective  value. 
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Duality  to  Solve  for  the  Penalty  Multiplier 


One  method  of  obtaining  a  lower  bound  on  the  true  optimal 
value  is  to  find  a  dual  feasible  point.  Using  the  gradient 
method  in  [1]  to  solve  the  Lagrangian  Dual  Problem  will  yield 
a  lower  bound  for  the  primal  problem  being  solved.  For  the 
example  problem,  the  Lagrangian  function  is: 

L  {x,  u)  -f(x)  u^gr^  (x) 

Then,  set  up  the  dual  problem: 

max  0  (u)  ,  u^O 


where: 


Q{u)  -  infimum  L(x,  u)  ,  X  €  01" 

Fortunately,  all  that  is  needed  for  the  primal  problem  is  any 
lower  bound  on  the  optimal  value.  Therefore,  the  dual  can  be 
initialized  with  u=(0,0),  thus  reducing  L(x,u)=f(x).  So,  to 
get  a  lower  bound  for  the  convex  primal  problem,  one  needs 
only  to  minimize  the  original  objective  function  without  any 
of  the  constraints. 

For  the  example  problem,  initializing  u=(0, 0,0,0),  will 
result  in  the  problem:  min  fix)  -  min  (x^-4 ) -+ (x2-4 ) ’ 

which  has  the  solution  x-(4,4)  with  f  ix)  -f'-O 
A  value  for  the  penalty  multiplier  can  now  be  established: 
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with  the  final  formulation  for  the  new  unconstrained 
minimization: 


min  (Xj^-4)  (X2-4)  ^+19  max  (x)  ,  0] 


NOTE:  Knowing  that  the  optimal  value  for  this  problem  is  8 

allows  calculation  of  a  smaller  multiplier  equal  to  9 .  In  the 
next  chapter,  both  of  these  values  will  be  used  to  solve  the 
problem. 
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CHAPTER  2 :  METHODOLOGY 


[5],  [6],  and  [9],  Ecker  and  Kupferschmid  develop  and 

implement  an  ellipsoid  algorithm  for  solving  nonlinear 
programming  problems.  For  a  given  n-dimensional  NLP,  the 
ellipsoid  algorithm  generates  a  sequence  of  shrinking 
ellipsoids.  Each  of  these  ellipsoids  contains  the  optimal 
point.  The  volume  of  the  ellipsoids  continues  to  reduce,  thus 
the  algorithm  converges  to  the  optimal  solution.  The  version 
of  the  algorithm  used  in  these  computations  was  EA3,  which  can 
be  found  in  [9]. 


The  Ellipsoid  Algorithm 

The  ellipsoid  algorithm  starts  with  an  initial  center 
point  and  initial  ellipsoid: 

Eq  -  [xe9l"  1  (x-x°)  ^1] 

where  x°  is  the  center  point  of  the  initial  ellipsoid  and  Qq 
is  a  symmetric  positive  definite  matrix.  This  initial 
ellipsoid  contains  the  optimal  solution  (guaranteed  if  this 
ellipsoid  contains  the  entire  feasible  set) .  The  algorithm 
generates  sequentially  smaller  ellipsoids  through  what  is 
termed  either  a  phase  1  or  phase  2  iteration  process. 

A  phase  1  iteration  occurs  if  the  center  point  of  the 
current  ellipsoid  is  infeasible.  A  cutting  hyperplane 
calculated  from  the  gradient  of  the  first  violated  constraint 
is  used  to  generate  the  next  ellipsoid.  A  phase  2  iteration 
occurs  if  the  center  point  of  the  current  ellipsoid  is 
feasible  (there  are  no  violated  constraints) .  For  this 
iteration,  the  cutting  hyperplane  is  calculated  using  the 
gradient  of  the  objective  function.  The  general  form  for  the 
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equation  of  the  cutting  hyperplane  (where  k  is  the  index  for 
the  current  iterate)  is: 

Hjf  -  [x  (x*)  -0] 

where:  i=  index  of  violated  constraint  for  a  phase  1  iterate 

m+l  (objective  fen  gradient)  for  a  phase  2  iterate 

In  both  cases,  the  current  ellipsoid  is  cut  in  half  and  a 
new  ellipsoid  is  generated  which  includes  the  half  of  the 
previous  ellipsoid  containing  the  optimal  point.  The  volume 
reduction  on  each  iteration  is  dependent  on  the  number  of 
variables  and  can  be  expressed  as  the  ration  of  volumes  of  two 
successive  ellipsoids: 


VoKE^,^) 

n-1 

n 

Vol(E^)  \ 

n-t-1 

(n2-i)i/2 

Some  example  values  for  the  amount  of  volume  reduction  on  each 
iterate  are: 

n=  #  of  variables  q^=  volume  reduction 


2  .77 

10  .9511 

100  .9950 

As  the  number  of  variables  in  the  problem  increases,  the 
volume  reduction  for  each  iterate  decreases.  In  general,  the 
fewer  the  variable  that  exist,  the  faster  the  convergence  is 
to  the  optimal  solution.  In  its  current  form,  the  ellipsoid 
algorithm,  EA3,  will  continue  to  run  until  either  the  optimal 
point  is  found  or  the  current  matrix,  Q|^  is  no  longer 
numerically  positive  definite. 
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Function  and  Gradient  Routines 


Specific  routines  exist  for  coding  nonlinear  programming 
problems  to  be  solved  by  a  specific  algorithm.  For  solving 
general  constrained  problems,  the  NLP  is  coded  into  two 
fortran  routines,  FCN  (for  function  calls)  and  GRAD  (for 
gradient  calls) .  These  two  routines  are  then  passed  as 
parameters  to  whichever  algorithm/program  is  being  used.  For 
the  example  problem  started  in  chapter  1,  the  function  (FCN) 
and  gradient  (GRAD)  routines  can  be  found  in  appendix  3.  Each 
of  the  problems  examined  were  likewise  coded. 

Some  additional  subroutines  are  required  to  solve  the 
exact  penalty  function  formulated  NLP.  By  it's  construction 
the  reformulated  problem  has  no  constraints,  yet  the 
constraint  components  are  included  in  the  objective  function. 
Also,  the  ellipsoid  algorithm  searches  for  a  violated 
constraint  function  before  turning  to  the  objective  function. 
In  the  case  of  the  reformulated  problem,  the  algorithm  need 
only  look  directly  to  the  gradient  of  the  objective  function 
(since  there  are  no  constraints  to  violate,  hence  no  gradient 
evaluation  for  a  violated  constraint) . 

Two  additional  subroutines  FCNX  and  GRADX  were  written  to 
take  advantage  of  the  existing  FCN  and  GRAD  routines,  yet 
perform  only  the  calculations  needed  for  the  exact  penalty 
function  formulation.  FCNX  prompts  for  the  original  number  of 
constraints  and  the  penalty  multiplier  to  be  used.  The 
routine  then  uses  the  FCN  routine  for  the  specific  problem  to 
construct  the  complete  unconstrained  exact  penalty  function. 
GRADX  uses  the  gradient  information  provided  in  the  problem 
GRAD  routine  and  constructs  the  gradient  for  the  unconstrained 
objective  function.  The  listings  of  FCNX  and  GRADX  can  be 
found  at  appendix  4 . 


10 


Error  Analysis 


Error  versus  effort  curves  are  generated  for  each  of  the 
various  problems  to  trace  the  accuracy  of  a  solution  against 
time.  The  horizontal  axis  is  the  effort  axis.  Effort  is 
measured  by  problem  state  central  processing  unit  (PSCPU) 
time.  The  only  time  counted  in  solving  the  problem  is  during 
the  steps  of  the  actual  algorithm  execution.  For  the  penalty 
function  problem,  since  no  constraints  exist  which  could  be 
violated,  then  there  is  no  need  to  calculate  function  values 
during  iteration.  The  gradient  of  the  objective  function  is 
required  on  every  iterate  so  only  those  calculations  are 
counted  for  time. 

The  vertical  axis  of  the  error  versus  effort  plot  is 
titled  "Log  Relative  Combined  Solution  Error."  The  combined 
error  for  each  iterate,  x*',  is  calculated  as  follows: 

e(x'')  -  I  fg  (x'") -fg  (x*)  I  +  I  I 

where  x*  is  the  optimal  point  and  the  are  the  Lagrange 

multipliers  at  optimality.  These  error  terms  are  then 
normalized  to  obtain  relative  combined  error  terms: 


E(x*^) 


e(x 

e(x°) 


The  logarithms  of  these  E(x‘')  error  terms  are  then  plotted 
against  the  PSCPU  times.  The  error  versus  effort  curves  for 
the  example  problem  is  at  appendix  2,  page  2. 
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CHAPTER  3!  RESULTS  AND  CONCLUSIONS 


Sample  Problems 

Each  of  the  problems  examined  were  first  solved  as 
constrained  NLP's  using  the  ellipsoid  algorithm,  EA3 .  Then, 
using  the  calculations  demonstrated  for  the  example  problem  in 
chapter  1,  a  penalty  multiplier  for  each  was  determined.  The 
problem  was  then  solved  as  an  unconstrained  exact  penalty 
function.  The  problems  examined  were: 

Table  I:  Problem  Summaries 

PROBLEM _ n=  #  VARIABLES _ m=  #  CONSTRAINTS _ REF 


Example 

2 

4 

N/A 

Colville  1 

5 

15 

[3] 

BBZ  3 

16 

13 

[2] 

Powell  19 

2 

20 

[4] 

Fukushima 

5 

11 

[8] 

For  the  example  problem,  the  strictly  feasible  point  used 
was  z=(l,l).  Two  lower  bounds  on  the  optimal  objective  value 
were  calculated  for  use  in  determining  a  penalty  multiplier. 
The  first  bound  was  based  on  the  unconstrained  minimization  of 
the  objective  function  (initial  iteration  of  the  gradient 
method  for  solving  the  Lagrangian  Dual  Problem) .  The  second 
was  based  on  knowledge  of  the  optimal  solution.  Two 
corresponding  penalty  multipliers  were  then  calculated.  The 
same  procedure  was  followed  to  calculate  penalty  multiplier 
values  for  the  remainder  of  the  sample  problems.  For  each  of 
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the  other  four  sample  problems,  the  procedure  is  only 
complicated  by  the  larger  number  of  variables  in  each  problem 
and  the  larger  number  of  constraint  equations.  The 
calculations  for  z  and  the  lower  bound  on  the  objective 
function  were  simplified  through  the  use  of  a  Fortran  program 
by  Mike  Kupferschmid,  AM46 ; FEASTEST.  This  program  examines  a 
given  point,  x,  and  determines  if  any  of  the  constraint 
equations  in  the  FCN  routine  are  violated.  The  program  also 
lists  the  constraints  in  order,  from  closest  to  being  violated 
to  farthest.  From  the  screen  output,  the  x  being  investigated 

is  clearly  strictly  feasible  or  not.  Also  the  a  value  can  be 

picked  right  from  the  rank  order  of  constraints.  Results  for 
these  preliminary  calculations  are  shown  here: 


Table  II:  Penalty  Multiplier  Calculations 


PROBLEM  Z=STRICT  FEAS  PT _ F  '  =LOW  BNP 


Example 

[1] 

Colville  1 

[.1  .1  .2 

BBZ3 

[.1  .6  0 

-3.7  -2. 
.7  5  0. 

Powell  19 

[-.5  -.5] 

Fukushima 

[012-1 

0  and  8 
-61.448 
30.62331 

-1.9999987 
-353 . 0495 


C=PEN  MULT 
19  and  9 
443.56 
41 

7 

11 


Optimal  Value  Results 

The  sample  problems  were  solved  using  the  ellipsoid 
algorithm  on  each  of  the  two  forms  of  the  problem  (constrained 
and  penalty  function) .  The  best  feasible  solution  point  was 
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then  found  for  each  solution  process  and  these  results  were 
then  tabulated  for  comparison: 


Table  III:  Optimal  Value  Results 


PROBLEM 
Example 
Colville  1 
BBZ  3 
Powell  19 
Fukushima 


CONSTRAINED  SOLUTION 
8.000000000000009 

-32.34867896572270 

30.62330955021475 

-1.414231127874633 

-43.99999999999998 


PENALTY  FCN  SOLUTION 

8.000000000000000 

-32.34867896572270 

30.62334552248174 

-1.414231127874633 

-43.99999999999998 


The  best  solution  is  shown  in  bold  for  each  of  the 
problems.  In  the  case  of  the  simple  example  problem,  the 
penalty  function  formulation  solved  the  NLP  to  the  exact 
solution  of  8.  Solution  of  the  constrained  problem  yielded  an 
optimal  value  within  lO"’®  of  the  exact  solution.  For  the 
problems  Colville  1,  Powell  19,  and  Fukushima,  both 
formulations  yielded  the  same  result  for  the  optimal  values. 

In  the  case  of  BBZ  3 ,  the  solution  of  the  exact  penalty 
function  formulated  problem  yielded  a  solution  within  10'^  of 
the  solution  to  the  constrained  problem.  These  results 
indicate  that  for  convex  problems,  the  solution  of  the 
reformulated  NLP  achieves  basically  equal  results  with  the 
constrained  solution. 


Time  to  Solution  Results 

The  other  part  of  the  comparison  of  the  methods  of 
formulation  is  to  see  if  the  reformulated  problem  can  be 
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solved  faster.  Intuitively,  the  solution  to  an  unconstrained 
problem  would  be  easier  or  faster  to  achieve  as  compared  to  a 
constrained  problem.  After  each  problem  was  solved  in  each  of 
the  two  formulations,  the  largest  time  difference  between 
iterations  was  calculated,  and  the  results  are  presented  here: 


Table  IV:  Solution  Time  Results 


PROBLEM 


Example 
Colville  1 
BBZ  3 
Powell  19 


LARGEST  TIME  DIFFERENCE  (%  FASTER) 
CONSTRAINED  PROBLEM _ PENALTY  FCN  SOLUTION 


14.6  % 
23.9  % 
54.2  % 
14.6  % 


Fukushima 


6.4  % 


These  results  show  that  in  four  of  the  five  problems,  the 
constrained  formulation  converged  to  a  solution  faster  than 
the  penalty  function  formulation.  Only  the  Fukushima  problem 
was  solved  faster  in  the  unconstrained  formulation.  The  error 
versus  effort  curves  in  appendix  2  support  these  results.  The 
plots  for  Colville  1  (pg  2-3) ,  BBZ  3  (pg  2-4) ,  and  Powell  19 
(pg  2-5)  all  show  steeper,  hence  faster,  convergence  for  the 
constrained  formulation.  Only  the  Fukushima  plot  (pg  2-6) 
shows  the  unconstrained  penalty  function  formulation 
converging  at  a  faster  rate. 

The  solution  time  results  run  counter  to  the  intuitive 
expectations  as  to  which  formulation  of  the  problem  should  be 
faster  to  solve.  An  explanation  for  these  results  comes  from 
the  manner  in  which  the  ellipsoid  algorithm  solves  a  problem. 

A  standard  NLP  has  "m”  constraint  equations  with  the  objective 
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function  labelled  as  the  ”in+l''st  equation.  Recall,  the 
ellipsoid  algorithm  evaluates  the  constraints  until  one  is 
found  which  is  violated.'  The  corresponding  gradient  for  that 
constraint  is  then  evaluated.  These  function  and  gradient 
evaluations  are  the  timed  segments  in  the  effort  portion  of 
the  plots.  The  number  of  overall  function  evaluations  can, 
therefore,  vary  on  each  iteration: 

1.  Best  Case:  the  first  constraint  examined  is  violated. 

The  result  is  that  only  a  single  function  and  its  gradient  are 
evaluated.  Total  function  evaluations:  2. 

2.  Worst  Case:  all  of  the  constraints  are  satisfied  (the 

current  point  is  feasible) .  So,  all  m  constraints  are 
evaluated.  The  objective  function  and  its  gradient  are 
evaluated.  Total  function  evaluations:  (m+l)+l=  m+2. 

3.  All  Others:  an  intermediate  constraint  is  violated,  so 
greater  than  two  evaluations  but  less  than  m+2  evaluations  are 
needed. 

The  exact  penalty  function  formulation  requires  no 
function  evaluations  from  FCN  to  solve  the  problem.  The 
cutting  hyperplane  is  determined  based  on  the  gradient  of  the 
only  function  present,  the  objective  function.  However,  since 
the  objective  function  includes  a  sum  of  all  of  the 
constraints  (the  p(x,c)  penalty  function),  the  gradient  of  all 
"m"  constraints  plus  the  original  objective  function  gradient 
are  calculated  on  each  iterate.  The  result  is  "m+1"  function 
evaluations  on  each  iterate.  No  variations  occur  from 
iteration  to  iteration.  Graphically,  the  number  of 
evaluations  per  iteration  can  be  shown  as: 


Bajoa 

e  for*  Constrained.  Forriu  lat,  ion _ 

1 - 

2 

Exact  Penalty  Fonnulation 

- T" 

1 

n*2 

Figure  2:  Function  Evaluation  Requirements 
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Conclusions 


The  exact  penalty  function  formulation  provides  an 
alternative  means  to  solve  convex  nonlinear  programming 
problems.  The  theorem  has  been  proved  (appendix  1)  to  always 
converge  to  a  solution  for  convex  problems,  achieving  the  same 
optimal  point  as  when  solving  the  standard  constrained  NLP, 
when  analytic  solutions  are  possible.  Differences  in  the 
solution  can  occur  due  to  numerical  methods  in  the  solution 
algorithm,  as  seen  in  the  optimal  solution  results  for  the 
problems  examined. 

Zangwill  presents  a  method  to  calculate  a  sufficient 
penalty  multiplier  to  guarantee  convergence  to  the  same 
optimal  point.  The  theorem  also  guarantees  convergence  for 
any  multiplier  greater  than  the  one  calculated.  However, 
this  penalty  multiplier  is  not  necessarily  the  minimum 
multiplier  which  will  guarantee  convergence.  For  the  simple 
example  problem  worked  throughout  this  project,  two  penalty 
multipliers  were  calculated  (9  and  19) .  The  error  versus 
effort  plot  showing  these  two  convergence  trajectories  is  at 
page  2  of  appendix  2.  The  multiplier  9  was  based  on  the  true 
objective  value  of  8  for  the  problem.  However,  lower  penalty 
multipliers  will  also  cause  the  solution  to  converge  until  a 
multiplier  of  3  is  chosen.  The  plot  for  convergence  using  a 
penalty  multiplier  of  3  is  also  shown  on  the  same  error  versus 
effort  plot.  The  plot  shows  that  the  solution  diverges  from 
the  true  optimal  value,  while  the  plots  for  multipliers  9  and 
19  basically  follow  the  same  path. 

In  four  out  of  five  test  problems,  the  reformulated 
unconstrained  problem  was  solved  at  a  slower  rate  than  the 
original  constrained  problem.  Only  in  the  Fukushima  problem 
was  faster  convergence  achieved  by  the  penalty  function 
formulation.  As  shown  in  figure  2,  the  penalty  function 
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formulation  requires  that  m+1  functions  be  evaluated  on  each 
iteration.  Meanwhile,  each  iteration  of  the  constrained 
problem  requires  evaluations  in  the  range;  [2,  m+2]. 
Therefore,  for  the  reformulated  problem  to  converge  faster, 
the  ellipsoid  algorithm  must  demonstrate  worst  case  behavior 
(every  iterate,  x*',  is  a  feasible  point)  .  This  worst  case 
behavior  is  demonstrated  in  the  Fukushima  problem  where 
virtually  every  iterate  is  feasible.  In  this  case,  the 
ellipsoid  algorithm  solves  the  reformulated  unconstrained 
problem  faster.  However,  this  type  of  behavior  in  the 
solution  process  is  more  the  exception  rather  than  the  norm. 
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APPENDIX  l:  PROOF  OF  2AMGWILL  THEOREM 


In  [11],  Zangwill  states  and  proves  his  penalty  function 
theorem  for  the  concave  nonlinear  programming  problem. 
Restated,  the  theorem  for  the  convex  nonlinear  problem  is  as 
follows: 

THEOREM:  let  X*  be  the  optimal  point  for  the  convex  nonlinear 

programming  problem: 

minimize  f (x) 

s.t.  gj(x)  ^  0  ,  i— l,2,...,m 
where  x  6  R"  and  f  and  the  gj '  s  are  convex 

and,  let  S°={x  j  g,.  (x)  <0,  i=l,2,...,m)  =  interior  of  the 

feasible  set 

and,  assume:  S°  is  nonempty,  and  S°=int(S) 

then,  the  convex  NLP  can  be  solved  by  a  single  unconstrained 
minimization 


PROOF: 

1.  define:  p (x,  c) -f  (x) (x)  ,  0] 

then,  minimizing  p(x,c)  yields  unconstrained  NLP 

2.  find  a  c  ,  so  that  3  min  [p(x,  c]] 

ie:  p(x*,  -n\injj[p(x,  c)  ]  where  x*=  optimal  point 

and: 

3.  by  assumption,  S°  is  nonempty  so  let:  z  6  S°, 
s.t.  gj(z)<0 

define:  a-max_;  ( z)  ]  <0 

4.  let  X*  be  the  optimal  point,  then  define: 

P-f(x*)-f(z)  and  c- 

a 
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5.  now,  given  any  point  w  ,  which  is  infeasible,  finding  a 
point  V  ,  feasible,  where:  p( v,  c)  <p(w,  c)  ,  is  sufficient  to 
prove  the  theorem. 

ie:  minp(x,c^  on  S~p{x* ,c) -f [x") 

-now,  for  any  w  infeasible,  show 
that  there  exists  a  v  feasible, 

s.t.:  p( V,  c)  <p(w,  c) 

-this  implies  that  the  minimizing 
point  for  p{x,c)  over  R"  must  be 
in  the  feasible  set. 


6.  define  v  as  follows  (see  figure  1.1):  z€S °  w$S 

which  implies  that  3v66(S')  ,  and  v  lies  on  the  line  segment 
joining  z  and  w 

now:  ve6  (S)  implies  3  set  A- [i  I  ( v) -0] 

and  Vi  C  A,  g^(v}  <0 

7.  define  a  function:  t  (x) -f  (x) +c  (x)  ,  convex 

so,  since:  g’f(v)-0,  VieA  —  t(v) -f  (v) +c  V.  g Av) -p{v,c)  ~f  (v) 

8.  at  w:  g^{w)^0,  VieA  (because  w  is  infeasible) 

® °  *  I^ieA  ^  ~Y^i€A  [  s^i  ( '  0  ]  ^Hr.i  ^  0  ] 

9.  — t(w) -f{w)  +C  V  gi(w)  ^f{w)  +C  V"  max  [g^  {x)  ,  0] 

or:  t  (w)  sp(w,  c) 
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*  The  proof  began  with  premise  of  showing  that  p{v,c)<p{w,^ 

but:  t(v) -p(v,  c)  =*t(v’)  <p(w,  c)  ^t(w) 

therefore,  to  complete  the  proof,  it  is  sufficient  to  show 
that:  t(v)  <  t(w) 

10.  first  show  t(z)  <  t(w): 

t(z)  -fiz)  -fiz)  + 5; g,  ( z) 

t  iz)  (z)  +  max.[gy(z)]  ,  because  gj(z)<0  for  all  i 

a 

t(z)  ^f(z)  +  ~^[^\~^max^[g^(z)  ] 

maxjg^(z)]  ^ 

t  (z)  ^f(x*)  -1 

t{z)<f{v)  ,  since  v  is  feasible,  z  is  strictly  feasible 

t  (z)  <  t  (v)  ,  since  t(v)=f(v) 

11.  now  show  t(v)<t(w): 

rewrite  v  as:  v-Xz+(l-X)w,  A.e(0,l)  ,  and  t(x)  convex 
=*t(v)^A,t(z) +  (!->.)  t{w) 

t(v)  <Xt(v)  +  (l-A.)  t(w)  ,  since  t(z)<t(v) 

(1-X)  t(v)  <  (1-1)  t{w) 

<t(w) 

12 .  then  X*  is  the  optimal  point  and  the  theorem  holds 
V  c^c 
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APPENDIX  2;  ERROR  VS  EFFORT  PLOTS 


On  the  following  pages  are  the  error  versus  effort  plots 
for  the  convex  nonlinear  programming  problems  examined  during 
this  study.  They  are  included  as  follows: 

1.  Simple  Example  Problem . 2-2 

2.  Colville  1 . 2-3 

3.  BBZ  3 . 2-4 

4.  Powell  19 . 2-5 

5.  Fukushima . 2-6 
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Log  Relative  Combined  Solution  Error 
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Total  ^SCPU  Time  (sec) 
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Log  Relative  Combined  Solution  Error 

r16.8Q  -14.40  -12.00  -9.CQ  -?.20  -4.80  -2.40  0.00  2.40 
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APPENDIX  3;  FUNCTION,  GRADIENT  ROUTINES 


C 

Code  by  Christopher  Fowler  RPI  Troy,  NY  12180 

C 

DOUBLE  PRECISION  FUNCTION  FCN(X,N,I) 

C  This  routine  computes  a  function  value  for  the  example 

C  problem  used  in  the  masters  project. 

C 

C  variable  meaning 

C  -  - 

C  I  index  of  function  whose  value  is  wanted 

C  N  number  of  variables  (=2) 

C  X  point  at  which  the  value  is  wanted 

C 

REAL*8  X(N) 

C 

C - 

Calculate  the  required  constraint  function  value 
C 

IF(I.EQ.l)  FCN=  X(1)-2.D0 
IF(I.EQ.2)  FCN=  X(2)-2.D0 
IF(I.EQ.3)  FCN=-X(1) 

IF(I.EQ.4)  FCN=-X(2) 

C 

Calculate  the  objective  function  value 
C 

IF(I.EQ.5)  FCN=(X(1)-4.D0) **2  +  (X (2) -4 . DO) **2 
C 

RETURN 

END 


C 

Code  by  Christopher  Fowler  RPI  Troy,  NY  12180 

C 


C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 


SUBROUTINE  GRAD (X , N , I , G) 

This  routine  computes  a  gradient  vector  for  the  example 
problem  used  in  the  masters  project. 

variable  meaning 


G  gradient  vector  returned 

I  index  of  function  whose  gradient  is  wanted 

J  index  on  the  variables 

N  number  of  variables  (=2) 

X  point  at  which  the  gradient  is  wanted 


REAL*8  X(N),G(N) 
C 

C - 


APP  3,  PAGE  1 


ono  o  o  o  ooo 


start  with  the  gradients  of  the  constraints 


IF(I.NE.l)  GOTO  2 

1  G(1)=1.D0 
G(2)=0.D0 
RETURN 

2  IF(I.NE.2)  GOTO  3 
G(1)=0.D0 
G(2)=1.D0 
RETURN 

3  IF(I.NE.3)  GOTO  4 
G(1)=-1.D0 
G(2)=  O.DO 
RETURN 

4  IF(I.NE.4)  GOTO  5 
G(l)=  O.DO 
G(2)=-1.D0 
RETURN 


now  calculate  the  gradient  of  the  objective  fen 

5  G(1)=2.D0*X(1)-8.D0 
G(2)=2.D0*X(2)--8.D0 
C 

RETURN 

END 

This  is  file  CONVEX 

It  contains  FCN  and  GRAD  for  an  example  nonlinear  programming 
problem  for  the  masters  project. 
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APPENDIX  4 


FCNX  AND  GRADX  SUBROUTINES 


I.  Function  Call  fFCNX^ : 


C 

Code 

C 

C 

C 

C 

C 

C 

C 

C 

C 

c 

c 


BLOCK  DATA 

COMMON  /BOUNDS/  XHIN , XLIN, NIN, MIIN, MEIN 
REAL*8  XHIN(50)/50*Z8181818181818181/ 

REAL*8  XLIN(50)/50*Z8181818181818181/ 

INTEGER*4  NIN/Z81818181/,MIIN/0/ , MEIN/0/ 

END 

by  Michael  Kupferschmid  and  Chris  Fowler 
SUBROUTINE  GET$IN 

This  routine  gets  the  penalty  mult  at  start  of  a  run. 
Variable  Meaning 


CBAR  penalty  multiplier 

FREAD  system  routine  for  free-format  input 

MEPF  number  of  original  constraints 

PROMPT  routine  prompts  for  input  from  the  terminal 

send  the  penalty  multiplier  from  common 
COMMON  /EPF/  CBAR, MEPF 
REAL* 8  CBAR 


prompt-  for  and  read  number  of  original  constraints 
CALL  PROMPT (' original  constraints= ' , 21) 

MEPF=0 

CALL  FREAD ( ' GUSER ' , • INTEGER*4 : ' , MEPF , & 1 ) 

prompt  for  and  read  the  penalty  multiplier 

1  CALL  PROMPT ( 'penalty  multiplier= ' , 20) 

CBAR=l.D+06 

CALL  FREAD ( 'GUSER' , 'REAL* 8; ' ,CBAR,&2) 

2  RETURN 
END 

by  Mike  Kupferschmid  and  Chris  Fowler 

DOUBLE  PRECISION  FUNCTION  FCNX(X,N,I) 

This  routine  computes  the  function  value  for  a 
penalty  function  problem 


C 
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C  Variable  Meaning 

C  - 

C  CBAR 

C  DMAXl 

C  FCN 

C  I 

C  II 

C  MEPF 

C  N 

C  X 

C 

REAL*8  X(N),FCN 

receive  the  penalty  mult  from  common 
COMMON  /EPF/  CBAR, MEPF 
REAL* 8  CBAR 


penalty  multiplier 

Fortran  fen  gives  larger  values  REAL*8's 
original  function  routine  for  problem 
unused;  index  of  the  pen  fen  objective (=1) 
index  of  original  fen  whose  value  is  wanted 
number  of  original  constraints 
number  of  variables  in  problem 
point  at  which  fen  values  are  wanted 


compute  the  value  of  the  exact  penalty  function 
objective  function 
FCNX=FCN ( X , N , MEPF+ 1 ) 

constraints 
DO  1  11=1, MEPF 

FCNX=FCNX+CBAR*DMAX1 ( 1 . DO , FCN ( X , N , I I ) ) 

1  CONTINUE 

RETURN 

END 


II.  Gradient  Call  fGRADX) ; 


C 

Code  by  Michael  Kupferschmid  and  Chris  Fowler 
C 


C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 

C 


SUBROUTINE  GRADX (X , N, I , G) 

This  routine  computes  the  gradient  vector  for  the 
penalty  function  problem 

Variable  Meaning 


CBAR  penalty  multiplier 

FCN  original  function  routine  for  the  problem 

G  gradient  vector  returned 

GC  gradient  vector  of  a  constraint  function 

GRAD  original  gradient  subroutine 

I  unused;index  of  pen  fen  objective (=1) 

II  index  of  orig  fen  for  which  gradient  needed 

J  index  on  the  variables 

MEPF  number  of  original  constraints 
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number  of  variables  in  problem 
Point  at  which  function  values  wanted 


N 
X 

REAL*8  X(N),G(N) 

REAL*8  GC(50),FCN 

receive  the  penalty  multiplier  from  common 
COMMON  /EPF/  CBAR,MEPF 
REAL* 8  CBAR 


compute  the  gradient  of  the  exact  pen  fen  objective 
CALL  GRAD(X,N,MEPF+1,G) 

constraints 
DO  1  11=1, MEPF 

check  whether  the  constraint  is  satisfied  at  X 
IF(FCN(X,N, II) .LE.O.DO)  GOTOl 

constraint  is  violated;  add  in  CBAR* (natural 
gradient) 

CALL  GRAD(X,N,II,GC) 

DO  2  J=1,N 
.  J,  =G(J)+CBAR*GC(J) 

2  CONTINUE 

1  CONT_NUE 
RE'T’URN 
END 
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